\section{Apéndice B: Código}

Incluímos en esta sección las funciones más relevantes del tp.

\subsection{Covarianza}

{\small 
\begin{verbatim}
Matriz covarianza(Album a, int ancho) {
    Iterador it;
    Matriz aux(ancho*ancho);
    double promedio[ancho*ancho];
    for(unsigned int i=0;i<ancho*ancho;i++) {
        promedio[i] = 0;
    }
    int cantidad=0;
    for ( it=a.begin() ; it != a.end(); it++ ) {
        for(unsigned int i=0;i<ancho;i++) {
            for(unsigned int j=0;j<ancho;j++) {
                promedio[i*ancho + j] += (double)(*it).second.pos(i,j);
                aux.setPos(cantidad, i*ancho+j, (*it).second.pos(i,j));
            }
        }
        cantidad++;
    }
    //cout << aux << endl;
    //hasta aca en promedio[] tenemos el prmediode las cosas
    for(unsigned int i=0;i<ancho*ancho;i++) {
        promedio[i] = (double)promedio[i]/cantidad;
    }
    for(unsigned int i=0;i<cantidad;i++) {
        for(unsigned int j=0;j<ancho*ancho;j++) {
            aux.setPos(i,j,aux(i,j) -  promedio[j]);
        }
    }
    //datos en la matriz a usar para hacer los calculos
    Matriz cov(ancho*ancho);
    for(unsigned int c=0;c<cantidad;c++) {
        for(unsigned int i=0;i<ancho*ancho;i++) {
            for(unsigned int j=i;j<ancho*ancho;j++) {
                cov.setPos(i,j, cov(i,j) + aux(c, i)*aux(c,j));
            }
        }
    }
    for(unsigned int i=0;i<ancho*ancho;i++) {
        for(unsigned int j=i+1;j<ancho*ancho;j++) {
            cov.setPos(j,i, cov(i,j));
        }
    }
    cov.producto((double)1/(cantidad-1));
    //cout << cov << endl;
    return cov;
}
\end{verbatim}}

\subsection{Autocalcular}
{\small
\begin{verbatim}
/* Obtiene los autovalores y autovectores de la matriz
 * utilizando el metodo QR. Deja en autovalores, los autovalores de la matriz
 * y devuelve la Matriz que representa la basde de autovectores corresponeinde
 * en el order de autovalores
 */
Matriz Matriz::autocalcular(double* autovalores, double delta=0.1) const{
    Matriz A = (*this);
    double deltaTemp = delta;
    //Iteramos mientras A no sea casi triangular superior
    Matriz autovectores(this->size());
    autovectores.setIdentidad();
    while(A.sumaInferior() > deltaTemp) {
        //SIENDO A=Q*R dejamos en Q la matriz ortogonal y en A la matriz triangular superior(R)
        Matriz Q = A.householder();
        autovectores = autovectores*Q;
        //Y APLICAMOS EL METODO A=R*Q
        A = A*Q;
    }
    //Dejamos en autovalores la diagonal de A, que se corresponde
    //con sus autovalores, en orden decreciente segun el absoluto de los mismos
    A.getDiagonal(autovalores);
    return autovectores;
}

/* Devuelve la suma de los elementos de los absolutos
 * de todos los elementos que esten por debajo de la diagonal de A */
double Matriz::sumaInferior() const {
    double val=0;
    for(unsigned int i=1;i<size();i++) {
        for(unsigned int j=0;j<i;j++) {
            val+= fabs((*this)(i,j));
        }
    }
    return val;
}
\end{verbatim}}



\subsection{Householder}
{\small
\begin{verbatim}
/* Aplica el metodo de householder a la matriz actual "A"
 * Deja en A la matriz triangular superior R y devuelve la
 * matriz ortogonal Q talque A=Q*R 
 */

Matriz Matriz::householder(){
    /*Para obtener el metodo lo que se hace es
      ir reflejando la matriz para eliminar cada columna
      utilizando las matrices de householder. Hacemos esto
      hasta obtener Qn*..*Q1*A = R. Esa es la R de la factorizacion,
      y a Q la obtenemos haciendo Q^t = Qn*..*Q1 */
    Matriz Q(this->size());
    Q.setIdentidad();
    for(unsigned int i=0;i<size()-1;i++) {
        Matriz Qi = this->getQi(i); //La matriz que pone en 0 la i-esima columna
        Q = Qi*Q;//Q_ = Qi*..*Q1
        (*this) = Qi*(*this); //THIS=R=Qi*..*Q1*A=Qi*R_{i-1}
        //elimino errores numericos
        //Dejo en 0 toda la i-esiama columna (salvo el primer elemento)
        for(unsigned int j=i+1;j<size();j++)
            this->setPos(j,i,0);
    }
    //Y por ultimo transponemos Q=Qn*..*Q1
    Q.transponer();
    return Q;
}
\end{verbatim}}

\subsection{getQi}
{\small
\begin{verbatim}
 /* Devuelve la i-esima Matriz de Householder, para poner en 0
 * los elementos de la i-esima columna debajo de la i-esima fila */
Matriz Matriz::getQi(int columna) const {
    /* El metodo consite en reflexiones para obtener una matriz
       ortogonal talque esa matriz por la Matriz A, deje a A con 0
       en la columna deseada. Para realizar esto lo que hacemos es,
       dado la columna y fila i, nos quedamos con el vector que 
       que va desde [aii,..,ain]. a este vector le restamos en su
       primer componente, la norma del vector. A este nuevo vector
       lo llamamos u=[aii/norm,ai(i+1),..,ain]^t..
       A partir de aqui para obtener la matriz de householder, lo que hacemos
       es realizar la operacion Qi = I-2/norm(u)^2*u*u^t, y extenedmos la matriz
       para dejarla de tamno de la orignial, con 1 en la diagonal de los nuevos lugares*/
    int tamano = size()-columna; //el tamano del vector u
    double vectorResultante[tamano]; //este va a ser u
    double normaVector=0;//esta va a ser la norma del subvector [aii,..,ain]
    bool seguir = true;
    for(int i=0;i<tamano;i++) {
        vectorResultante[i] = (*this)(columna+i,columna);
        normaVector += vectorResultante[i]*vectorResultante[i];
    }
    //esto es FUNDAMENTAL, aca es donde se hacia las cosas con 0 que quedaban horrible
    //cuando la norma vector = al cuadrado del primer elemento del vector(ie, no hay q reflejar nada)
    //entoces en este caso se devuelvce la matriz Identidad! Veo dsps como resolver esto del 0.00001
    if((normaVector)-(vectorResultante[0]*vectorResultante[0]) < 0.0001)  {
        Matriz ret(size());
        ret.setIdentidad();
        return ret;
    }
    //aca vectorResultante = [aii,..,ain]
    normaVector = sqrt(normaVector); //normaVecotr = ||[aii,..,ain]_2]||
    vectorResultante[0] = vectorResultante[0] - normaVector; //vectorResultante = u
    //ahora calculamos la norma2 al cuadrado
    //de este vector y la dejamos en normaVector
    normaVector = 0;
    for(int i=0;i<tamano;i++) {
        normaVector += vectorResultante[i]*vectorResultante[i];
    }
    //Creamos Qi que es la matriz a devolver, esta tiene
    //entre de 0..columna-1:1 en la diagonal, 0 furera
    //entre columna..size() la submatriz I-2/norm(u)^2*u*u^t
    Matriz Qi(size());
    Qi.setIdentidad();//tenmos la I, ahora falta restar
    for(int i=0;i<tamano;i++) {
        for(int j=0;j<tamano;j++) {
            double restar;
            restar = vectorResultante[i]*vectorResultante[j];//esto es u*u^t en la pos adecuada
            restar = restar*2/normaVector;
            Qi.setPos(columna+i,columna+j, Qi(columna+i,columna+j)-restar); //aca queda para cada punto
        }
    }
    //y listo
    return Qi;
}
\end{verbatim}}

\subsection{mejorDistancia}
{\small
\begin{verbatim}
 string mejorDistancia(Album a, double * test, int cantCoordenadas, string * grupos, int cantGrupos) {
    //PARA CADA GRUPO CALCULAMOS SU DISTANCIA MAHALANOBI Y NOS QUEDAMOS CON LA DE MENRO DISTANCIA
    double distancia[cantCoordenadas];
    for(int i=0;i<cantGrupos;i++) {
        pair<Iterador, Iterador> group = a.equal_range(grupos[i]);
        double promedioGrupo[cantCoordenadas];
        Matriz cov = covarianzaTC(group.first, group.second, cantCoordenadas, promedioGrupo);
        cout << cov << endl;
        cov.inversa();
        for(int j=0;j<cantCoordenadas;j++) {
            test[j] = test[j] - promedioGrupo[j];
        }
        double values[cantCoordenadas];
        //cout << test[0] << "," << test[1] << " <<----" << endl;
        for(int j=0;j<cantCoordenadas;j++) {
            double val = 0;
            for(int k=0;k<cantCoordenadas;k++) {
                val += test[k]*cov(k,j);
            }
            values[j] = val;
        }
        //cout << values[0] << " , " << values[1] << "  lala <<" << endl;
        double valFinal = 0;
        for(int j=0;j<cantCoordenadas;j++) {
            valFinal += values[j]*test[j];
        }
        distancia[i] = sqrt(valFinal);
        cout << "FINAL(" << grupos[i] << "): " << distancia[i] << endl;
    }
    int indice = 0;
    for(int i=0;i<cantGrupos;i++) {
        if(distancia[i] < distancia[indice]) {
            indice = i;
        }
    }
    return grupos[indice];
}
\end{verbatim}}
